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When a ship suffers underwater damage, there is a rapid influx of water, followed by a 
period of slower progressive flooding. This results in flooding of compartments whose hull 
boundaries, but not interior bulkheads, are still intact. An existing computer model uses the 
FORTRAN computer language and formatted input files to model progressive flooding. 
This thesis uses MATLAB computer language and SIMULINK graphical user interface to 
provide a modular, expandable progressive flooding design tool. 
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I. INTRODUCTION 



A. MOTIVATION 

Accurate modeling of progressive flooding is essential in assessing the survivability 
features of a given ship. A realistic flooding and counter flooding scenario is a very 
complex dynamic system including complicated vessel motions, time-varying external 
excitation, changing hull mass and hydrodynamic properties, progressive flooding 
through damage openings, internal sloshing of flooded water, and de-watering control 
actions. Although, due to the complexity of the problem, the traditional method of 
approach is through large amplitude static analyses, the need for a more accurate analysis 
has been recognized. In order to bridge this gap, recently there has been a large amount of 
research focusing on dynamic stability prediction taking into consideration water-on-deck 
and large amplitude time-dependent motions and excitation [Ref. 1], Equally important, 
however, is the development of a design methodology for incorporating elements of 
man/machine interface into the problem; in other words effects of decision-making and 
control strategies. This presents an additional level of difficulty since combinations of 
discrete events and continuous differential equations describe such systems. The long 
term goal that is initiated with this thesis is to study the development of an improved 
damage stability analysis and design tool, which will address the effects of control 
decisions and vessel dynamics in a progressive flooding scenario. 
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B. BACKGROUND 



When a ship suffers underwater damage, there is a rapid influx of flooding water, 
usually followed by a period of slower progressive flooding. Progressive flooding is 
defined as flooding which extends beyond the compartments that are open directly to the 
sea. This results in flooding of compartments whose hull boundaries, but not interior 
bulkheads, are still intact. Existing computer models at the Naval Postgraduate School 
utilize the FORTRAN computer language and formatted input files to model progressive 
flooding [Ref. 2]. Any change in model design conditions, flooding damage, or damage 
control abilities requires a completely new FORTRAN input file. The purpose of this 
thesis is to include the capabilities of the original computer progressive flooding model in 
a more modular and expandable format, using the user-friendly MATLAB and 
SIMULINK programs. 

C. DEFINING THE MODEL 

1. Hole Model 

Holes can be modeled as a short-tube orifice as shown in figure 1 with a diameter 
(hole size) much larger than its tube length (ship’s skin thickness). The flow discharge 
coefficient for this case is 0.816 [Ref.3: pp.42-43]. 

2. Hull Model 

The required hull characteristics are determined using the same table of offsets 
used by Dawson (Hull.dat) [Ref. 2: pp. 127-143] and MATLAB programs arvol . m and 
drafts . m to calculate the locations of center of buoyancy and center of floatation, 
waterplane area, and waterplane moments used in the program equations. Both Matlab 
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programs are included in the appendix. For the FFG-7 hull model used here the 
following values were calculated: 

Table 1: Hull Model Characteristics 
Center of Buoyancy 



Xbar 


209.9 ft aft of FP 


Ybar 


0.0 ft from CL 


Zbar 

Center of Floatation 


10.52 ft from baseline 


Xbar 


227.4 ft aft of FP 


Ybar 


0.0 ft from CL 


Waterplane area 


14030 ft 2 


Longitudinal Moment of Inertia.... 


1.5883 x 10 9 ft 4 


Transverse Moment of Inertia 

3. Equilibrium Model 


1.8148 x 10 6 ft 4 



The added weight method [Ref.4: pp.76-79] is used within the program to 
determine the equilibrium condition (displacement, heel, and trim) of the hull model for 
any state of flooding. This method assumes hull characteristics are constant and is 
accurate for “small angles”. The added weight method treats flooding water as added 
weight and treats the hull as if it were intact. The effect of the added weight is divided 
into a weight, which causes the ship to immerse evenly (parallel sinkage) and moments 
which cause the ship to list and trim about its center of floatation. The small angle 
assumption assumes that the ship lists and trims about a fixed center of floatation. In 
reality, since the hull is not cylindrical, the location of the center of floatation will change 
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with changing list and trim. This method allows much simplified calculation of hull 
characteristics compared with the lost buoyancy method. Since the weight added is liquid 
water which is free to move as the ship moves (and thus does not contribute to pitch and 
roll resistance), a free surface correction is added to the longitudinal and transverse 
moments of inertia to account for the reduced resistance of the hull to pitch and roll. 

4. Compartment Model 

The ship compartments are modeled as rectangular boxes. Extension to more 
realistic shapes is one of the recommendations for further expansions to the capabilities 
of the program. 

D. SIMULINK 

The SIMULINK program provides a graphical representation of mathematical 
relationships and equations to enable the user to determine the response of system 
odels. Graphical symbols and connections are used instead of numbers and 
mathematical symbols, allowing more intuitive use and understanding of complex 
mathematical relationships between system components. The SIMULINK program also 
allows easy modification of input parameters and import and export of data from 
MATLAB programs. Lastly, SIMULINK allows the user to observe system or individual 
component behavior during the simulation through the use of output graphs. 

E. OBJECTIVES 

The objective of this thesis is to develop a progressive flooding model to utilize the 
superior user interface and data exchange capabilities of the SIMULINK program. The 
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program must be able to accurately calculate flooding rates and the subsequent response 
of the ship. The program must also be user-friendly, easily modified and have the 
potential for further expansion by subsequent users. To achieve these goals the following 
programs were written by the author: 

1. Programs 

a. Pflood.m 

Pflood.m is the main program created by the author. This program calculates the 
flooding and ship response caused by a hole of user-chosen size and location in the 
ship’s hull. The flooding rate, flooded volume, heel and trim are calculated and 
updated each second. In addition, the user can specify the number, time of 
activation and size of pumps used for dewatering and the duration of their use. 

b. Arvol.m 

Arvol.m is a MATLAB program written by the author to determine the center of 
buoyancy of the hull model defined by the Hull.dat data file using Simpson’s rule 
integration techniques. 

c. Drafts.m 

Drafts. m is a MATLAB program written by the author to determine the waterplane 
area, center of floatation, and second moment of inertia of the hull waterplane for 
a given draft. 

2. Validation 

Program performance was validated by comparing the program output to expected 
ship behavior for a given flooding condition. 
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II. SIMULINK PROGRAM 



The schematic representation of Pflood.m is shown in Figure 1 . For illustrative purposes, 
Pflood.m is shown in its “expanded” form. By leaving Pflood.m in expanded form, any 
function or variable can be observed or changed by the user. For ordinary use, the various 
subroutines would be grouped into subsystem blocks with only the user-defined variables 
and the outputs represented by discrete function blocks. Placing the program in 
subsystem block form has the added advantage of greatly enhancing program readability 
and user-friendliness. Specific details of the individual block diagram elements are 
included in Chapter HE. 
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Figure 1: Pflood.m Simulink Program 
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III. PROGRAM AND SUBROUTINES 



A. FLOW RATE THROUGH HOLE(S) 

The flow through the damage, either primary or secondary, can be modeled as 
turbulent flow through an orifice of negligible length, or a hole. Flow rate through a hole 
is a function of the hole size (cross sectional area), hole shape, and pressure difference 
across the hole. The hole size and shape are related to the discharge coefficient, Cj. A Ca 
of 0.816 is used in this application in order to remain consistent with Dawson’s previous 
computer flooding model. This number is also consistent with theoretical predictions of 
flow through orifices as explained in the previous chapter [Ref. 3: pp. 42-43]. The 
pressure difference is determined by the difference in height of the liquid on either side of 
the hole. The equation is 

Q = CdxAx ^/2gx(/r2-/n) 

where Q is the volumetric flow rate, A is the cross sectional area, Cd is the discharge 
coefficient, g is the acceleration due to gravity, and h 2 -hi is the liquid height difference 
across the hole. This height difference is a function of the level of the flooding water 
inside the flooded tanks, as well as the draft and trim of the ship. The latter quantities 
will have to be calculated using the added weight method as explained in the previous 
chapter. The Simulink representation of the above equation is shown in Figure 2. 

The program allows continuous adjustment of the values of h 1 and h2 to account for 
increased immersion of the exterior hole and increased compartment water height. 
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Figure 2: SIMULINK Hole Flow Symbology 
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B. NET FLOW INTO A COMPARTMENT 

Flooding water removed from the compartment by pump activation or by flow 
into adjacent compartments is subtracted from the water flow rate into the compartment. 
Two way flow due to water removal from one of two connected flooded compartments is 
also possible, and is accounted for by taking the square root of absolute value of the 
liquid height difference and then multiplying the result by the sign of the height 
difference. This method avoids program failure due to trying to take the square root of a 
negative number. This program also continuously updates the values of hi and h2 to 
account for increased immersion of the exterior hole and increased compartment water 
height. 

C. HEIGHT OF WATER (HOW) 

Once the net volumetric flow rate of water is known, integrating the flow rate 
results in the net volume of water added to the compartment. This volume, divided by the 
horizontal area of the compartment, and divided again by the permeability of the 
compartment, gives the height of water in the compartment. A saturation function ensures 
that the height of water in the compartment is always between zero and the maximum 
height of the compartment. Figure 3 shows the Simulink height of water function. 

D. MOMENTS DUE TO ADDED WATER 

The added weight method treats the added weight as a combination of a weight 
which causes the ship to sink evenly (parallel sinkage) and transverse and longitudinal 
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Figure 3: SEMULINK Height of Water Terminology 
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moments which cause the ship to heel and trim respectively. Multiplying the height of 
water (HOW) by the product of the weight density of the liquid ( pg ), the compartment 
permeability (per), and the compartment horizontal area (A) determines the weight added 
to the compartment. Multiplying the added weight by the moment arm (Zu-ans/iong) from 
the ship’s center of floatation to the center of gravity of the added weight gives the 
moment about the center of floatation in the longitudinal and transverse directions. In 
equation form: 

Moment — ( pgAxHOW x per) X Ztruns / limp 
The Simulink representation of the above equation is shown in Figure 4. 

E. KG (VERTICAL POSITION OF CENTER OF GRAVITY) 

KG (Vertical Position of Center of Gravity) is determined by multiplying the 
original weight of the ship by the original KG, adding the product of the weight of the 
added water and the height of the added water center of gravity above the keel, and 
dividing the sum by the new ship weight. The equation for KG is: 

KGnew = (KG old x Wold + KG add x Wadd) Wnew 
The Simulink representation is shown in Figure 5. 

F. VERTICAL CENTER OF BUOYANCY (KB) 

The ship’s new vertical center of buoyancy (KB) is determined by multiplying the 
original KB by the original displacement, adding the product of the added weight and the 
original draft plus one half of the parallel sinkage, then dividing the product by the new 
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Figure 4: SEMULINK Moment Terminology 
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Figure 5: SIMULINK KG Terminology 
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weight. In equation form: 

KBnew = ( KBorig x Woric + KB add x Wadd) a- Wnew 
This is based on the assumption that the hull is effectively wall-sided, and that the vertical 
center of buoyancy of the newly submerged portion is located at one half the added draft. 
Both of these assumptions are consistent with the small angle assumption that was 
mentioned in the previous chapter. If the small angle assumption is no longer true, then 
the true position of the center of buoyancy will have to be calculated. This can be done by 
employing numerical integrations and utilizing the provided table of offsets, or by 
interfacing the program with a standard external ship hydrostatics calculations program, 
such as SHCP [Ref. 2: pp. 24-25]. Such an interface is relatively straight forward since 
the computational environment of Matlab/Simulink already allows incorporation of 
external ASCII files. Figure 6 shows the Simulink representation of the KB equation. 

G. HEIGHT OF THE METACENTER ABOVE THE CENTER OF BUOYANCY 
(BM) 

The height of the metacenter above the center of buoyancy (BM) is determined by 
the second moment of the waterplane (about the waterplane centerline or longitudinal 
center of floatation) divided by the hull displacement. For the added weight method, a 
free surface correction is made by subtracting the second moment of the flooded 
compartment area from the waterplane moment. Typically, free surface corrections are 
done for liquids when calculating the final metacentric height. In such a case, the free 
surface correction appears as a reduction in the metacentric height or, equivalently, an 
increase in the virtual center of gravity. The same effect can be observed by modifying 
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the second moment of the waterplane area in the calculation of the metacentric radius, 
BM. This is because the added weight in this case is seawater. With this modification in 
BM, the added weight method produces identical results with the lost buoyancy 
calculations (reference a naval architecture book). The BM Simulink terminology is 
shown in Figure 7. 

H. RIGHTING ARM (GM) 

The righting arm (GM) is determined by GM = KB + BM - KG. This equation is 
graphically represented by the triple summing point shown in Figure 8. 

I. HEEL AND TRIM ANGLE 

The angle of heel is determined by 

sin (j> = Volume x Ztrans ( Displacement x GM ) 
as shown in Figure 9. The total trim is determined by dividing the total trimming moment 
by the moment to change trim (MCT). MCT is equal to the total weight times the 
longitudinal GM, divided by the ship length at the waterline. The tangent of the trim 
angle is equal to the total trim divided by the ship’s length at the waterline. If the large 
angle assumption needs to be utilized, the moment to change trim will have to be 
imported as an external array. This is not expected to cause significant difficulties since 
MCT is readily available from the hydrostatic diagrams of a given ship. The Simulink 
trim angle terminology is shown in Figure 10. 
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Figure 7: SIMULINK BM Terminology 
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Figure 8: SIMULINK GM Terminology 
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Figure 9: SIMULINK Heel Angle Terminology 
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Figure 10: SIMULINKTrim Angle Terminology 



J. ADDITIONAL SUBMERGENCE 

The additional submergence of the hole due to list is determined by multiplying the 
sine of the heel angle by the transverse distance from the transverse center of floatation to 
the hole. Similarly, the increased hole depth due to trim is obtained by multiplying the 
sine of the trim angle by the longitudinal distance to the center of floatation. The new 
hole depth is the sum of the original hole depth and the depth increases due to heel and 
trim. In algebraic form, these relations are shown by 

h 2 new = h2oLD + ( Ztrans .Me xsin <j>) + ( Zlonc .Me X sin ( trim )) 
and in graphical form they are represented by the diagram of Figure 1 1 . 



K. PROGRAM USAGE 

To use Pflood.m, the user needs only select “simulation” from the SIMULINK 
menu, then select “start”. Modifications to the simulation parameters or the program are 
easily accomplished by “double-clicking” on the symbol of interest. Additionally, 
SIMULINK allows the user to easily add output displays in order to observe the behavior 
of any desired variable during the simulation. 
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Figure 11: SEMULINK Hole Submergence Terminology 
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IV. RESULTS 



Results for a test case consisting of a compartment with a 4 square foot hull 
penetration (compartment 1) and a 4 square inch watertight bulkhead penetration 
(compartment 2). Pumps were activated at 10 minutes (600 seconds) which completely 
dewatered compartment 1 in order to test the ability of the program to account for 
backflow through holes. As shown in figure 12, the four square foot hull penetration 
completely floods compartment 1 essentially instantly. Activating the pumps 
dewaters compartment 1, instantly removing the water entering through the hull 
penetration. Similar results could have been achieved by changing the characteristics of 
the hull penetration to simulate hull repairs. 

Figure 1 3 illustrates how the program accurately predicts backflow into 
compartment 1 from compartment 2, which had no pumps assigned. A saturation function 
ensures that backflow from compartment 2 will stop when compartment 2 water height 
drops to the level of the bulkhead penetration. 

The predicted trim response of the ship is shown in Figure 14. The ship initially 
trims as a result of the rapid flooding of compartment 1 , then gradually increases its trim 
angle as compartment 2 floods. The instant dewatering of compartment 1 then causes the 
trim angle to decrease, followed by further trim angle reduction as compartment 2 is 
dewatered. Without altering the pump or hole configuration, the ship will assume a steady 
trim angle corresponding to compartment 2 flooded up to the height of the bulkhead 
penetration. 
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Figure 12: Height of Water in Compartment 1 
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Figure 1 3: Height of Water in Compartment 2 
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Figure 14: Ship Trim Response 



V. DISCUSSION, CONCLUSIONS, AND RECOMMENDATIONS 



A. DISCUSSION 

This thesis has successfully developed the foundation for a modular and 
expandable progressive flooding computer model. The computer model accurately 
predicts the quasi-static ship response to progressive flooding and is easily modified by 
the user to explore ship response to any desired damage situations. 



B. CONCLUSIONS 

Based on the above results, the following conclusions are drawn: 

1 . Pflood.m is an effective design tool for modeling and studying the effects of 
progressive flooding. 

2. Pflood.m is easily modified to model virtually any damage scenario. 

3. Pflood.m can be expanded to include different ship progressive flooding and 
equilibrium models. 



C. RECOMMENDATIONS 

The following items are recommended to further expand Pflood.m and enhance 
its utility as a design tool: 

1. Incorporate a hull characteristics look-up table into the Pflood.m program 
to expand the ship response modeling ability beyond “small angles” of 
heel and trim and non-rectangular compartments. 

2. Explore the effect of pump configuration on dewatering effectiveness and 
ship response. 
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APPENDIX A: HULL.DAT 



-.98 .203 34.011 
-.98 1.0097 37.035 
-.98 2.09 40.059 
-.98 2.30593 40.636 
-.98 2.525 41.213 
-.49 .204 24.201 
-.49 1.0381 28.1305 
-.49 2.2 32.06 
-.49 3.60937 36.059 
-.49 5.199 40.058 
-.49 5.26214 40.2255 
-.49 5.325 40.393 
.196 .207 10.687 
.196 .4439 12.3765 
.196 .741 14.066 
.196 1.2739 17.065 
.196 1.881 20.064 
.196 2.9 24.061 
.196 4.249 28.058 
.196 5.8592 32.0575 
.196 7.529 36.057 
.196 8.2218 37 . 6 ^ 
.196 8.921 39.323 
.49 .205 5.841 
.49 .319 6.453 
.49 .431 7.065 
.49 .72482 8.5655 
.49 .981 10.066 
.49 1.3109 12.0655 
.49 1.651 14.065 
.49 1.847 15.0645 
.49 2.031 16.064 
.49 2.941 20.062 
.49 4.1 24.06 
.49 5.6094 28.058 
.49 7.338 32.056 
.49 8.8498 35.4695 
.49 10.346 38.883 
.98 .201 .831 
.98 .32937 1.193 
.98 .558 1.555 
.98 .84406 2.309 
.98 1.036 3.063 
.98 1.2508 4.064 
.98 1.441 5.065 
.98 1.621 6.065 
.98 1.811 7.065 
.98 2.10725 8.565 
.98 2.381 10.065 
.98 2.7609 12.0645 
.98 3.166 14.064 
.98 3.38543 15.0635 
.98 3.621 16.063 
.98 4.70985 20.0605 
.98 6.099 24.058 
.98 7.78322 28.0565 
.98 9.628 32.055 
.98 11.0459 35.1215 
.98 12.466 38.1880 

1.569 .443 .031 

1.569 .73779 .288 

1.569 .98 .545 
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1.569 


1.18262 .7975 


1.569 


1.365 1.05 


1.569 


1.667 1.553 


1.569 


1.918 2.056 


1.569 


2.30982 3.0585 


1.569 


2.63 4.051 


1.569 


2.911 5.062 


1.569 


3.171 6.063 


1.569 


3.41112 7.0635 


1.569 


3.641 8.064 


1.569 


4.081 10.064 


1.569 


4.521 12.064 


1.569 


4.79912 13.238 


1.569 


5.088 14.412 


1.569 


5.88294 17.2355 


1.569 


6.82 20 . 059 


1.569 


8.389 24.057 


1.569 


10.198 28.055 


1.569 


12.118 32.055 


1.569 


14.043 36.055 


1.569 


14.3628 36.724 


1.569 


14.682 37.393 


1.961 


.587 .031 


1.961 


.93975 .286 


1.961 


1.257 .541 


1.961 


1.53453 .7925 


1.961 


1.791 1.044 


1.961 


2.223 1.547 


1.961 


2.585 2.05 


1.961 


3.14751 3.054 


1.961 


3.569 4.058 


1.961 


3.91984 5.0595 


1.961 


4.22 6.061 


1.961 


4.491 7.062 


1.961 


4.741 8.063 


1.961 


5.221 10.063 


1.961 


5.701 12.063 


1.961 


6.0018 13.2375 


1.961 


6.31 14.412 


1.961 


7.17812 17.235 


1.961 


8.199 20.058 


1.961 


9.82922 24.0565 


1.961 


11.684 28.055 


1.961 


1 - 3.638 32.055 


1.961 


15.598 36.055 


1.961 


16.007 36.893 


2.353 


.593 .031 


2.353 


1.09723 .2845 


2.353 


1.528 .538 


2.353 


1.89568 .7895 


2.353 


2.217 1.041 


2.353 


2.776 11.544 


2.353 


3.233 2.047 


2.353 


3 .936 3.051 


2.353 


4.466 4.055 


2.353 


4.8888 5.0575 


2.353 


5.254 6.06 


2.353 


5.56 7.061 
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629 


19.1 


20 


.196 


3 . 


779 


19.168 


20 


. 196 


5 . 


634 


19.343 


20 


.196 


6 . 


956 


19 . 519 


20 


.196 


8 . 


664 


19.806 


20 


.196 


10 


.1 2 


> 0.09 


20 


.196 


11 


.796 20.478 


20 


. 196 


13 


.204 20.862 


20 


.196 


13 


.976 21.074 


20 


. 196 


14 


.746 21.285 


20 


.196 


15 


. 05 


22.114 


20 


. 196 


14 


.34 


22 . 943 


20 


.196 


15 


.78 


24.397 


20 


.196 


16 


.17 


25.851 


20 


.196 


16 


.569 27.4 


20 


.196 


16 


.947 28.951 


20 


.294 


.584 21.361 


20 


.294 


4 . 


413 


21.536 


20 


.294 


6 . 


005 


21.711 
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20.294 8.225 22.049 

20.294 10.009 22.386 

20.294 13.064 23.1275 

20.294 15.414 23.869 

20.294 16.149 26.406 

20.294 16.884 28.943 

20.392 .588 23.67 

20.392 4.155 23.851 

20.392 5.975 24.031 

20.392 8.307 24.3715 

20.392 10.009 24.712 

20.392 12.188 25.2275 

20.392 14.012 25.743 

20.392 14.953 26.0365 

20.392 15.762 26.33 

20.392 16.284 27.637 

20.392 16.805 28.943 

20.49 .534 26.011 

20.49 3.117 26.086 

20.49 4.004 26.161 

20.49 7.627 26.599 

20.49 10.009 27.037 

20.49 12.125 27.535 

20.49 13.929 28.033 

20.49 15.086 28.4075 

20.49 16.131 28.782 

20.49 16.506 28.8625 

20.49 16.88 28.943 
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OP OP OP dP 



APPENDIX B. ARVOL.M 



arvol.m 

by Russell wright 

calculates station areas and center o' bouvancy 
for T =16.9 ft 



load hulll.dat 
a=zeros (41,3) ; 
zz=16 . 9 
for i=l : 41 



beg=23* <i-l)+l; 

endd=23 *i ; 
a ( i , 1) =hulll (beg) ; 
z=hulll (beg : endd, 3 ) ; 
if z ( 1 ) > zz 

a ( i , 3 ) =0 ; a ( i , 2 ) =1 ; 
else 

% y=hulll (beg : endd, 2 ) ; 

zp=hulll (beg : endd, 3 ) ; 

Z=zp(l) : ( 16 . 9-zp ( 1 ) )/22:16.9; 
yy=hulll (beg : endd, 2 ) ; 
y=sDline ( zd , yy , z ) ; 
h=z ( 2 ) -z ( 1 ) ; 
hhh(i) =h; 

four = 0 . 0 ; f r=0 ; 
two=0 . 0 ; tO=0 ; 

for j= 1:11 

four=f our+y (2* j ) ; 
f r=f r+z ( 2* j ) *y(2* j ) ; 

end 

for k= 1:10 

two= two+y (2*k+l ) ; 
to=to+y ( 2 *k+l ) *z(2*k+l) ; 

end 

% calculate station area using simpson's rule: 
a(i,2)=2My(l) +4 # four+2*two+y ( 23 ) ) *h/ 3 ; 

% calculate station centroid (z bar) using simpson's rule 
a ( i , 3 ) =2 * (y ( 1 ) *z ( 1 ) +4 * fr+2* to+y ( 23 ) *z { 23 ) ) *h/ ( 3 *a ( i , 2 ) ) ; 

end 

end 

aa=zeros (101 , 3 ) ; 

hh= ( (a ( 37 , 1) -a(3, 1) ) /100) ; 

aa ( : ,l) = (a(3,l) :hh:a(37,l) ) ' ; 

aa( : , 2 ) =spline ( a ( 3 : 37 , 1 ) , a (3 : 37 , 2) , aa ( : , 1) ) ; 

aa ( : , 3 ) =spline (a ( 3 : 37 , 1 ) , a ( 3 : 37 , 3 ) , aa ( : , 1 ) ) ; 

f ore=0 ; too=0 ; f re=0 ; tw=0 ; f o=0 ; tu=0 ; 
for i=l : 50 
fore=fore+aa (2*i, 2) ; 
fre=fre+aa (2*i,2) *aa ( 2 * i , 3 ) ; 
f o=f o+aa ( 2* i , 2 ) *aa (2 * i , 1 ) ; 
end 

for i=l:49 

too = too+aa (2*i+l,2) ; 
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tw= 

tu= 

end 

vol = 
zbar 
xbar 



tw+aa( 2 *i+l, 2 ) *aa( 2 *-+i 3 ) . 
tu+aa(2*i-i-l,2) •>aa( 2 *i+i'i) • 



*(aa(l,2)+4’fore*2*coo+aa(101,2) ) 

{ aa (x , 2 ) *aa ( X , 3 ) +4 * f ire+2 
= (aa ( 1 , 2 ) *aa (1,1) +4 



*cw+aa(101,2) *aa(101,3) ) *hh (3-voi) 
fo+2*tu+aa (101,2) *aa ( 101 , 1 ) ) *hh/ (3 ‘volt 
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APPENDIX C. DRAFTS.M 



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

%% 

% arafts2.m 

% by russell wright 

% calculates wp points for a given draft 

% 



load hulll.dat 
ind= 0 ; 

yincl = zeros (41,3) ; 
zz= 16 . 9 ; 



zincl = z z 
zinc2=zz 
for i=l : 41 
i 

beg=23* (i-1) +1; 
endd=23 *i ; 

wpla ( i , 1 ) =hulll (beg, 1 ) ; 
z =hu 1 1 1 ( beg : endd , 3 ) ; 
y=hulll (beg : endd, 2 ) ; 

wpla ( i , 2 ) =interpl ( z , y , zz , ' spline ' ) ; 
if z ( 1 ) > zz 
wpla ( i , 2 ) =0 ; 
yincl ( i , 2 ) =0 ; 
yincl ( i , 3 ) =0 ; 
wpla ( i , 1) =0 ; 
end 

yincl ( i , 1) =hulli (beg, 1 ) ; 
yincl (i, 2) =interpl (z,y, zincl, ' spline' ) ; 
if yincl(i,2) <0.0 
yincl (i, 2) =0.0; 
end 

yincl ( i , 3 ) =zz ; 

if wpla ( i , 2 ) ==0 
yincl ( i , 2 ) =0 ; 
yincl ( i , 3 ) =0 ; 
end 

end 

yincl=yinci (3:37, : ) ; 



xx=vincl ( : , 1 ) ; 

x=yincl(l, 1) : (yincl (35, 1) -yincl (1,1) ) /34 : yincl (35,1) ; 
yy=yincl { : , 2 ) ; 
y= spline (xx, yy , x) ; 



h=x ( 2 ) — x ( 1 ) ; 
hhh ( i ) =h; 

four = 0.0; fr=0 ; 
two=0 . 0 ; to=0 ; 



for j= 1:17 

four=four+y (2 X j ) ; 
fr=fr^x{2 # j ) *y (2* j ) ; 
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end 



for k=l:l6 

two= two+y (2 *k+l ) ; 
to=to+y (2*k+l) *x(2*k+l) ; 

end 

% calculate wp 1/2 area using simpson's rule: 
a2=2‘* (y (1) +4*four- l -2* two+y (35) ) *h/3 ; 

% calculate wd centroid (x bar) using simpson's rule: 
xbar=2’ (y (1) *x( 1) +4* fr+2*to+y (35) *x ( 3 5 ) ) *h/ ( 3 *a2 ) 

% calculate long, moment o' inertia about xbar 



f r=0 ; 
to = 0 ; 

for j= 1:17 

fr=fr+ (x(2*j)~2)*y(2*j) ; 

end 

for k=l:16 

to=to+y (2*k+l) * (x(2*k+l) A 2) ? 

end 

io=2* (y ( 1 ) *x ( 1 ) *x(l)+4*fr-+‘2* to-+*y (35) *(x(35)~2))*h/3; 
isubl = io-a2*xbar~2 

% calculate trans moment o' inertia 



f r=0 ; 
to=0 ; 

for j= 1:17 

fr=fr+(y(2*j)~3) ; 

end 

for k=l:16 

to=to+(y(2*k+l)~3) ; 

end 

isubt = 2 * ( (y(l) /N 3)+4*fr+2*to+(y(35)''3) ) *h/9; 
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